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Type I migration in optically thick accretion discs 
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ABSTRACT 

We study the torque acting on a planet embedded in an optically thick accretion 
disc, using global two-dimensional hydrodynamic simulations. The temperature of an 
optically thick accretion disc is determined by the energy balance between the viscous 
heating and the radiative cooling. The radiative cooling rate depends on the opacity 
of the disc. The opacity is expressed as a function of the temperature. We find the 
disc is divided into three regions that have different temperature distributions. The 
slope of the entropy distribution becomes steep in the inner region of the disc with the 
high temperature and the outer region of the disc with the low temperature, while it 
becomes shallow in the middle region with the intermediate temperature. Planets in 
the inner and outer regions move outward owing to the large positive corotation torque 
exerted on the planet by an adiabatic disc, on the other hand, a planet in the middle 
region moves inward toward the central star. Planets are expected to accumulate at 
the boundary between the inner and middle regions of the adiabatic disc. The positive 
corotation torque decreases with an increase in the viscosity of the disc. We find that 
the positive corotation torque acting on the planet in the inner region becomes too 
small to cancel the negative Lindblad torque when we include the large viscosity, 
which destroys the enhancement of the density in the horseshoe orbit of the planet. 
This leads to the inward migration of the planet in the inner region of the disc. A 
planet with 5 Earth masses in the inner region can move outward in a disc with the 
surface density of 100 g/cm 2 at 1 AU when the accretion rate of a disc is smaller than 
2 x 10- 8 A</ o /yr. 

Key words: hydrodynamics - radiative transfer - planets and satellites: formation - 
planetary systems: gravitational interactions - planetary systems: accretion disc. 



1 INTRODUCTION 



Planets are thought to be formed in a protoplanetary disc 
around a young star. A planet growing by the accretion of 
planetesimals exchanges the angular momentum with disc 
gas a nd moves in the disc |Goldreich fc Tremaind 1 19791 . 
1 19801 ). This is called the Type I migration of a planet. The 
Type I migration is caused by the torques acting on a planet 
by a disc. The torque is composed of the Lindblad torque 
and the corotation torque. The Lindblad torque is due to 
the two spiral den sity waves excited by the planet in the disc 
|Wardlll986l . ll997l ), while the corotation torque is exerted by 



the gas in the horse shoe r egion of the planet dWardl 



Baruteau fc Massetl 120081 



iPaardekooper fc Mellemal 



1991 



2008 



Paardekooper fc Papaloizoull2008l ). The sum of the Lindblad 



torque and the corotation torque determines the total torque 
acting on the planet. 

The Lindblad torque is exerted by the disc gas at the 
Lindblad resonances, which are located inside and outside 
the orbit of the planet. The angular momentum of the planet 
is increased and decreased by the inner and outer Lindblad 
torques, respectively. The magnitude of the negative outer 
Lindblad torque is a little larger than that of the p ositive 
inner Lindblad torque in a disc (|Wardlll98rj . 1 19971 ). The 
Lindblad torque becomes negative, leading to the inward 
migration of a planet. 

The corotation torque exerted by the disc gas in 
the horseshoe orbit of a planet depends on the vorten- 
sity distribution and the entropy distribution of the disc 
gas. The vortensity related corotation torque is small 
compared with the Lindblad torque (|Masset fc Casolil 
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IPaardekooper et alj l201ll ). on the other hand, the 
entropy related corotation torque can become larger 
than the magnitude of the negative Lindblad torque. 
IPaardekooper fc Mellemal (|2006l ) showed that the corotation 
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torque becomes positive in a disc with the high opacity, 
reducing the migration velocity of a planet. This process 
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itive corotation torque cancels the negative Lindblad torque 
when the entropy distribution has steep negative slope. It 
was also shown that a planet in an adiabatic disc with the 
steeper negative slope starts to migrate even outward. 

The corotation torque significantly decreases in an adi- 
abatic disc after a few libration time of a planet because the 
entropy distribution becomes uniform within the horsesho e 
region |Masset fc Casolil |2010| ; iPaardekooper et al.l l201ll ). 
Without the sustained corotation torque, only the Lindblad 
torque is eventually exerted on a planet in an adiabatic disc, 
leading to the inward migration of a plane t with Earth mass 
in 10 5 yrs (|Wardll 19971 ; iTanaka et alj|2002h . This is a serious 
problem in the core accretion model of planet formation. 
The lifetime of a disc, 10 6-7 yrs, is much longer than the 
timescale of the Type I migration, making the survival of 
planets in a disc difficult. However, recent radial velocity 
surveys of extrasolar planets show that a significant frac- 
tion of solar-type s tars may harbor clos e -in super-Earths 
(iMavor et al.l 120091 ; lLo Curto efaLI |2010| ; ISegransan et al.l 
uOllf h Population synthesis models have great difficulties to 
reproduce the observed semimajor axis distribution of ex- 
trasolar planets once the Type I migration is included. The 
reduction of the m i gration velocity of a p lanet is required 
|Alibert et alj|2005l ; iMordasini et al.ll20od ). 

The viscosity o f a disc prevents the s at uration of the 
corota tion torque (jMasset fc Casolil |2010T |. iKlev fc Cridal 
(|2008| ) investigated the planet-disc interactions in an op- 
tically thick accretion disc, taking into account the viscous 
heating and the radiative cooling. They found the disc gas 
reaches a steady state, of which the temperature distribu- 
tion has the very steep slope in a region where ice exists 
in a condensed form. The corotation torque becomes always 
positive in this region and the direction of the plan et migra- 
tion becomes outward. lHasegawa fc Pudritj (|201ll ) compre- 
hensively examined the various mechanisms to halt planet 
migration in a disc around a young star. They showed that 
the planet migration is halted at the ice line, inside of which 
all the ice is evaporated and solid particles are composed 
of rocks and metals. The slope of the temperature distribu- 
tion becomes steep in the inner region. The positive corota- 
tion torque becomes large enough to make the total torque 
positive. On the other hand, the temperature distribution 
beyond the ice line is too shallow to cancel the negative 
Lindblad torque by the corotation torque. A planet inside 
and outside of the ice line is expected to move outward and 
inward, respectively. They showed that the ice line is the 
location to stop the planet migration in an optically thick 
accretion disc. 

iMuto fc Inutsukal (|2009h showed that the torque is de- 
pendent on the magnitude of the viscosity. The viscosity 
modifies the dens ity of gas near a planet and decr eases 
the torque density. IPaardekooper fc Papaloizoul (|2009aT l also 
showed that the corotation torque becomes very small when 
they consider a disc with high viscosity. In this study, we sys- 
tematically examine the total torque acting on a planet by 



an optically thick accretion disc. It is not clear if the positive 
corotation torque can always exceed the negative Lindblad 
torque when we include viscosity in a disc. We make global 
two-dimensional hydrodynamic simulations to study the ef- 
fect of the dissipation processes such as the viscosity and 
the radiation on the total torque. 

This paper is organized as follows. In the section 2, we 
describe the basic equations and a disc model we presume 
in this study. The dissipative terms due to the viscosity and 
the radiation are included in the basic equations. The tem- 
perature of the disc is determined by the energy balance 
between the viscous heating and the radiative cooling. The 
rate of the radiative cooling is dependent on the opacity of 
the disc, which changes around the ice line. In the section 3, 
we show the results of the two-dimensional hydrodynamic 
simulations. We show the corotation torque decreases with 
an increase in the viscosity. The total torque acting on the 
planet depends on the radiative cooling rate and the vis- 
cous heating rate. We further derive the analytical formula 
to relate the viscosity and the opacity when the total torque 
exerted on the planet becomes zero. We summarize the re- 
sults in the section 4. 



2 BASIC EQUATIONS AND DISC MODEL 
2.1 Basic Equations 

A planet excites density waves in a disc and changes the 
density distribution of the disc. We examine the torque ex- 
erted on a planet by an optically thick accretion disc. A 
planet with 5 Earth masses rotates around a solar mass star 
in a fixed circular orbit. The position vector of the planet 
from the star is denoted by r p . The problem is limited to 
a two-dimensional flow, where all physical quantities (e.g., 
the surface density) depend on r and 9, where r is the dis- 
tance from the star and 9 is the angle between the x-axis 
and the position vector. Governing equations for the gas are 
the mass conservation, the Navier-Stokes equations, and the 
energy equation with dissipative terms due to the viscosity 
and the radiation. 

We use a cylindrical coordinate where the star is located 
at the center of the coordinate. The mass of the planet is 
much smaller than that of the star and we neglect the indi- 
rect term. The mass conservation equation and the Navier- 
Stokes equations read 



1 9 , v, s 1 d ^ s 



0. 



r dr 



r dr 



Ev r v e E <9$ 

+ J 9 

r r o9 



(1) 



(2) 



(3) 



where £ is the gas surface density, p is the vertically in- 
tegrated pressure; v r and Vg are the radial and tangential 
velocities of the gas; $ is the gravitational potential. We con- 
sider a less massive disc and neglect the self-gravity of the 
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gas. The gravitational potential of the star and the planet 
is given by 

GMq GM p 



4> 



y/ r 2 + r p — 2rr p cosifj + e 2 H p 



(4) 



where Mq is the mass of the sun, M p is the mass of the 
planet, ip is the angle between r and r p , and H p is the scale 
height of the disc at the location of the planet. The scale 
height is given by 



V2c p 



(5) 



where c p and f2 p are, respectively, the isothermal sound 



velocity and the Keplerian angular velocity at 



The 



smoothing length parameter, e, is introduced to include 
the effect of the scale height of a disc. It is noted that 
the small so ftening parameter leads to the strong corota- 
tion t orque jBaruteau fe Masset] [2008: Paa rdekooper et al.l 
l20ld ). lYamada fc Inaba I ( |201ll ) examined the torque acting 



on a planet with several Earth masses in an optically thin 
disc and found the torques calculated by the simulations 
with e = 0.3 agree with that obtained by a linear analy- 
sis. We adopt e — 0.3 in this study. The last terms in the 
Navier-Stokes equations describe the radial and azimuthal 
components of the viscous forces: 



19. 1 Oa r e 

-— (rorr) H 57- 

r Or r 06 



fr = 
and 



Q'ee 
r 



Or 



1 

+ - 

r 



(6) 



(7) 



where oy r , a r $, and age are the components of the stress 
tensor and are, respectively, given by 



2T,u 



Or 



and 



2Eu (dv e 



+ V r 



(8) 
(9) 

(10) 



The gas is assumed to be ideal and the equation of state is 
given by 

Efe B T 



Limn 



(11) 



where ks is the Boltzmann constant, T is the mid- plane 
temperature, mu is the mass of a hydrogen atom, and is 
the mean molecular weight of the gas. We set up fi =2.34. 
The energy equation of the gas disc reads 

9(Ee) 19. _ 1 9 , 



dp v e dp 
Vr d-r + Vd6 +W 



(12) 



where e is the specific energy of the gas given with the pres 
sure by 

P 



(7-1)2" 



(13) 



In Eq. (|13|) . 7 is the ratio of the specific heats at constant 
pressure and volume. We use 7 = 4/3 for a two-dimensional 



disc as given by I Li et all (|200Cf ) . The energy of the gas is 
increased by the viscous heating term. The viscous heating 
term, W, is given by 



w = 



a 2^L+a [r ( Ve \ I 1 ) 
rr dr r \ dr \ r J r 09 ) 

1 dvg v r 



+o"se 



/l dvg Vr_\ 
\r 09 r ) 



(14) 



On the other hand, the energy of the gas is decreased by the 
radiative cooling term of the disc, Q. The radiative cooling 
term, Q, is given by 



Q = 2CTSBr off , 



(15) 



where ass is the Stefan-Boltzmann constant and T e s is an 
effective temperature of a disc. The radiative flux in a radial 
direction is neglected because it is much smaller than Q by 
a factor of H p /r p . The effective tempera ture is related to 
the mid-plane temperature of the disc as (|Hubenvlll990h 



with 



3r 

Teff = — 



1 

4? 



(16) 



(17) 



where r is the optical depth of a disc. The optical depth of 
a disc is defined by 



T = -Cgr^oE, 



(18) 



where k,q is the be st fit function of the opacity 
jLin fc Papaloizou|[l985T ) and £ gr i s the s o-called grain con- 
tent factor introduced by iMizun 3 (|l980l ). Micron size dust 
particles are the main source of the opacity. Collisions be- 
tween particles might produce a large numb er of small dust 
parti cles, increasing the opacity of a disc (|Birnstiel et all 
l2010t ). The opacity might decrease due to capture of small 
dust particles by large particles. It is valuable to study the 
torque acting on a planet in a disc with various values of £ gr . 
We treat £ gr as a parameter of the simulations and consider 
the discs with 0.1 ^ £ gr ^ 100. 

The opacity strongly de pends on the temperature of a 



Ihe opacity strongly dej 
disc (|Lin fc Papaloizo jfl985h 



K0 = 



5 x 10~' 5 T 



2 x 10 _4 T 2 



(210K < T < 2000K), 
(170K < T < 210K), 
(T < 170K). 



(19) 



Equation (|19|l was derived from the available opacity data, 
which incorporates small particles found in interstellar 
medium. Dust particles are composed of rocks and metals in 
a region of disc with T ^ 210 K, while ice is added to form 
dust particles in a region of disc with T < 170K. In the tran- 
sition region between the cold and the hot region, the ice is 
evaporating. We divide a disc into three regions: the region 1 
with T > 210 K, the region 2 with 170 K < T < 210 K, and 
the region 3 with T < 170 K. 

For the later convenience, we write all the quantities 
in a non-dimensional form using the unit length ro= 1 AU, 
the unit mass Mq, and the unit time J7q 1 where Qo is the 
Kepler frequency at 1 AU. The normalized quantities are 
denoted by a tilde (e.g., f p ). 
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2.2 Disc Model 



iKitamura et all (120021 ) obtained images of protoplanetary 
discs around T Tauri stars in Taurus using the thermal emis- 
sion of dust. They found that the surface density distribu- 
tions of protoplanetary discs are described by the power-law 
distribution with the power-law index of 0-1. The surface 
densities of the discs at 100 AU are in a range between 0.1 
and 10 g/cm 2 . In this study we adopt a disc model, of which 
the initial condition is given by the power-law distribution: 



E ini = Eof" 



(20) 



where Eo is the surface density at ro and is set to be 

100 g/cm 2 . 

IKitamura et all l|2002h also found disc radii expand with 
time, suggesting that most of the gas migrates toward in- 
ward a central star by transporting the angular momentum 
outward . We adopt an accretion disc model developed by 
IPringld jl98ll ). in which the accretion rate of gas is con- 
stant throughout a disc. The viscosity v is expressed by a 
power-law function of the distance from the central star to 
satisfy the constant accretion rate of a disc as 



(21) 



where vo is the kinematic viscosity at ro and we set vq = 
4.2 x 10 15 cm 2 /s and £ v is a viscous strength factor. The value 
of i/o is approxim ately equivalent to . 05 in terms of the al- 
pha coefficient of IShakura fe Sunvaevl(|l973t ). The accretion 
rate of a disc with the kinematic viscosity, vo, and the surface 
density given by Eq. becomes 6.3 x 10~ 8 A/ o /yr. Obser- 
vations of discs sug gest that the accretion r ates of discs are 
~ 1 x 10" 8 M Q /yr jHartmann et al.lll998l ; IKitamura et all 
I2002T 1. We treat £ v as another parameter in addition to £ gr to 
consider discs with various accretion rates of 0.1 ^ £ v ^ 1. 
Two parameters, £ gr and £ v , determine the structure of a 
disc. 

We derive the temperat ure distribution of a n optically 
thick accretion disc following lKlev fc Cridal (|2008). The tem- 
perature distribution of the disc is determined by the bal- 
ance between the viscous heating term, W , and the radiative 
cooling term, Q. Assuming that the viscous heating is due 
to the Keplerian motion, we have W given by 

(22) 



W 



-Ei/fr 



Applying Eqs. (|16I) - (I19I ) to Q with the assumption of t > 
1, we obtain the temperature distribution of the disc from 
W = Q as 



210 H- 

210 (-!-) 

\ ri2 I 

170 (-?-) 



-(a+3)/3 

-(a+3)/ll 

-(a+3)/2 



(210K < T < 2000K), 
(170K sC T < 210K), 
(T < 170K), 



(23) 



where ri2 and r23 are the distances from the star to the 
boundary of the regions 1 and 2 and that of the regions 2 and 
3, respectively. For later convenience, the power-law indexes 
of the temperature distribution in the i-region are denoted 
by -ft (e.g., ft = (a + 3)/3). 

Figure [1] shows the surface density distribution and the 
temperature distribution of the disc with a — 1.0. It is seen 
from this figure that the temperature distribution has three 
different slopes. The boundary between the regions 1 and 2 
and that between the regions 2 and 3 are located at 1.4 AU 



and 2.5 AU, respectively. The temperature distribution is 
nearly flat in the region 2 (ft = 4/11). The power-law in- 
dexes of the temperature distributions are given by ft =4/3 
in the region 1 and ft = 2.0 in the region 3. 

The larger viscosity and/or opacity increase the tem- 
perature of the disc. The boundary positions move outward 
with increases in the viscosity and/or opacity. The location 
of the boundary is expressed by the two parameters, £ gr and 
£ v , and the boundary position between the regions 1 and 2 
is given by 



ri2 



1.4 



1/4 



(24) 



The boundary position between the regions 2 and 3 is given 
by T23 = l-8f 12 in the case of a = 1.0. This boun dary po- 
sitions were derived bv lHaseeawa fc Pudritj (|201lf ) as well. 
It is noted that the boundary position also depends on the 
surface density. We plot the boundary positions of the re- 
gions 1 and 2, fi2, as a function of £ gr and £ v in Figf2] The 
solid, dashed, dot-dashed, and dotted curves correspond to 
fi2 = 2.4, 1.4, 1.2, and 0.8, respectively. It is found that the 
boundary position between the regions 1 and 2 ranges from 
0.5 AU to 2.5 AU. The filled circles represent the parameter 
sets (£ g r,£v) used in our numerical simulations. 



3 NUMERICAL METHOD AND SIMULATION 
RESULTS 

A planet generates density waves in a disc inside and out- 
side of the orbit of the planet. The inner and outer density 
waves exert the positive and negative torques on the planet, 
respectively. The gravity of an outer density wave is a little 
stronger than that of an inner density wave because an outer 
density wave is closer to the planet due to the pressure gradi- 
ent o f a disc, leading to the negative Lindblad torque (|Wardl 
Il997h . Another torque acts on the pla net by a gas element 
in th e horseshoe orbit of the planet (iBaruteau fc Massetl 
120081 ; iPaardekooper fc Papaloizoul |200S? T When a gas ele- 
ment in the horseshoe orbit approaches a planet, the angu- 
lar momentum is exchanged between the gas element and 
the planet. This is called the corotation torque. 

The entropy of the gas is expressed as S — p/E 7 oc r A , 
where A = (7 — l)a — ft. The corotation torque is depen- 
dent on the entropy distribution of a disc. It was shown 
that the corotation torque consists of linear and non-linear 
corotation torque. The non-linear corotation torque comes 
from the outgoi ng boundaries of the horseshoe region, that 
is the separatrix (|Masset fc Ca soli 2010; Paa rdekooper et al.l 
|2010|) . In a non-barotropic disc, the vortensity is changed af- 
ter the horseshoe U-turn. The enhanced corotation torque 
is not due to the adiabatic compression, but comes from a 
singular streamline at the separatrix. Since the change of the 
vortensity is proportional to the radial entropy gradient of 
a disc, the large corotatio n torque is induced in a d isc with 
the large magnitude of A. lYamada fc Inaba! (|201 lh showed 
that the positive corotation torque is comparable with the 
negative Lindblad torque when A = —0.4 in an adiabatic 
disc. The total torque acting on a planet is determined by 
the sum of the negative Lindblad torque and the corotation 
torque. 

The disc we consider in the section 2.2 has the surface 
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density distribution with a = 1.0 and the temperature dis- 
tribution with pi = 4/3 in the region 1, /3a = 4/11 in the 
region 2, and p-j, = 2.0 in the region 3. The power- law indexes 
of the entropy distribution are calculated as A = —1.0 in the 
region 1, A = —3.0 x 10~ 2 in the regio n 2, and A = —1.7 in 
the region 3. lYamada fc Inabal (|201ll ) found that the total 
torque acting on a planet in an adiabatic disc becomes posi- 
tive when A < —0.4. We apply this to the optically thick disc 
and predict that a planet migrates outward in the regions 1 
and 3 while it moves inward toward the star in the region 2. 
Planets in the regions 1 and 2 are expected to move toward 
the boundary of the regions 1 and 2. This might help form- 
ing larger planets. However, it is not clear whether a planet 
still moves outward in the region 1 even if dissipative pro- 
cesses are included in a disc. We make a number of numerical 
simulations to examine the torque acting on a planet in a 
disc with dissipative processes due to the viscosity and the 
radiation. 



3.1 Numerical Method 

We use two-dimensional equidistant grids in r and 6 with 
a resolution of (N r , Ng) = (640, 3072). The inner and 
outer radii of the disc are given by r m i D /r p =0.4 and 
r m ax/ fp = 2.0, respectively. Th e damping boundary condi- 
tions l|de Val-Borro et alj|2006h . where all components are 
relaxed towards their initial state, are used in order to reduce 
wave reflection from these boundaries. All the quantities in 
the inner and outer boundaries are always fixed to be the 
initial values. 

We develop a two-dimensional global hydrodynamic 
computer program with the gravit ational forces of a star an d 
a planet and a dissipative term (|Yamada fc Inabal l201ll ). 
The basic equations are solved simultaneously using the fi- 
nite volume method with an operator splitting procedure. 
The source terms, which include the gravity, the viscous 
force, and the radiation, are computed with a second order 
Runge-Kutta scheme. The advection terms are calculated 
with a second or der MUSCL-Hancock scheme and an exact 
Riemann solver (|Torolll999l ; llnaba et al1l2005h . 

The angular momentum of the disc gas is transferred 
to the planet. The transfer rate of the angular momentum 
from the disc gas at r to the planet is given by 



9$ 

E w rde - 



(25) 



By integrating the torque density over the radial distance, 
we obtain the total torque acting on the planet from the 
disc: 



r = 



IVdr. 



(26) 



The torque and the torque density are normalized by To and 
r /r p , where T is (M p /A/ Q ) 2 (r p /H p ) 2 E p r P ^. 

3.2 Simulation Results 

3.2.1 The torque density acting on a planet in an 
optically thick accretion disc 

We consider two adiabatic discs that have the tempera- 
ture distributions with the power-law indexes 4/3 and 4/11, 



which correspond to that of the region 1 and the region 2 of 
the optically thick disc, respectively. The power-law index 
of the temperature distribution is expressed as — p. Both 
discs have the initial surface density with a = 1.0. Here- 
after we call the adiabatic disc with (a,P) = (1.0,4/3) and 
that with (a,P) = (1.0,4/11) the disc A and the disc B, 
respectively. The power-law index of the entropy distribu- 
tion of the disc A is A = —1.0 and that of the disc B is 
A = -3.0 x 10" 2 . lYamada fc Inaba I (|201ll) found that the 
corotation torque increases with a decrease in the power in- 
dex of the entropy distribution in an adiabatic disc. Figure[3] 
shows the radial distribution of the torque density exerted 
on the planet in the adiabatic discs at t = 25t p , where t p is 
the rotational period of the planet. The planet is located at 
1 AU in the disc A and 2 AU in the disc B. The large corota- 
tion torque is exerted on the planet in the disc A due to the 
large entropy gradient, making the total torque positive. On 
the other hand, the corotation torque in the disc B is too 
weak to cancel the negative Lindblad torque. The normal- 
ized total torques become 3.0 and —2.4 in the disc A and 
the disc B, respectively. 

The half width of the horseshoe region, Sx s , is estimated 
as (IMasset et al. 20061; Paardekooper fc Papaloizoul l2009bl ; 



iPaardekooper et al.ll201ol ) 



2.1 x 10" 




(27) 



where Me is the Earth mass and H p is the scale height of 
the disc at f = f D : 



— - (ft 



(l-/3)/2 



(28) 



iMasset et all (|2006l ) showed that the half width of the horse- 

1/2 

shoe region of a planet scales with M p as long as the 
flow around a planet remains linear. Figure U shows the half 
width of the horseshoe region given by the numerical sim- 
ulations (the filled circles) together with the half width of 
the horseshoe region given by Eq. ()27|) (the cur ve) . Both half 
width s agree well except the large planet mass. IMasset et al.l 
(|2006f ) found the planet mass when the flow linearity breaks 
in a 2D simulation. The non-linearity of the flow around a 
planet becomes clear when the mass of the planet becomes 
larger than 12A4e(# p /0.057- p ) 3 . The numerical results also 
start to diverge from the half width given by Ea. ()27p when 
the mass of the planet becomes larger than 10 Earth masses. 

The ratios of the magnitude of the corotation torque to 
that of the Lindblad torque are 1.6 and 0.6 in the disc A and 
the disc B, respectively. Nearly the same Lindblad torques 
act on the planet in both discs. The corotation torque de- 
termines the magnitude and sign of the total torque exerted 
on the planet. 

We perform the simulation to find the torque density 
exerted on the planet in the optically thick accretion disc 
with (£ gr ,£v) = (10,0.1). The temperature distribution is 
the same as that of Fig. [1] The boundary position between 
the regions 1 and 2 is 1.4 AU and that between the regions 2 
and 3 is 2.5 AU. Figure[5]shows the radial distributions of the 
torque density acting on the planets in the region 1 (panel 
(a)) and in the region 2 (panel (b)) at t = 25i p . The planets 
are located at 1 AU in the region 1 and 2 AU in the region 2. 
The radial distribution of the torque density in the region 1 



6 K. Yamada and S. Inaba 



is similar to that in the disc A. The torque density has a local 
maximum and a local minimum at |1 — r/r p \ = 3.3 x 10 -2 . 
The distances from the planet to the local maximum and 
minimum can be approximated by the scale height of the 
disc in the region 1, indicating that this torque corresponds 
to the Lindblad torque. 

The flow pattern within the horseshoe region is modified 
by the viscosity and the vortensity is no longer conserved 
during the U-turn. The vortensity is radially transferred by 
the viscosity. We compare the vortensity distributions of the 
inviscid disc, that is the adiabatic disc, and the accretion disc 
with non- vanishing viscosity. Figure |S] shows the contour of 
the vortensity within the horseshoe region of the planet. 
It is seen from Fig. [6] that the location of the minimum 
vortensity moves inward toward the central star when the 
viscosity is included. The viscosity plays an important role 
to diffuse the vortensity. In a disc with the low viscosity, the 
corotation torque is dominated by the non-linear corotation 
torque. The non-linear corotation torque decreases with the 
increasing viscosity. 

The radiation restores the entropy distribution of a disc 
and is i mportant to prevent the saturation of the corotation 
torque. iMasset fc Casolil (|20ld ) derived the analytical ex- 
pression of the corotation torque acting on a planet by a 
disc with viscosity and thermal diffusivity. In this study, the 
cooling due to the radiation is compensated by the viscous 
heating. The radiation only works to maintain the thermal 
structure of a disc. Our numerical results are compared with 
the analytical expressio n of the corotation torque given by 
IMasset fe Casolil (|2010l ') in Fig. \7\ We remove the contribu- 
tion of thermal diffusivity in the analytical expression. Our 
results agree well with the analytical result. In this study, the 
radiation affects the corotation torque exerted on a planet 
through the disc temperature profile. However, it is impor- 
tant to note that radiation significantly influ e nces t he coro- 
tation torque as pointed by IMasset fc Casolil (|2010h . 

The timescale for the gas to spread across the horseshoe 
width by the viscosity is defined by 

r visc = M . (29) 



We substitute Eqs. (J^TJ and (gT} into Eq. ((SSI to estimate 
the viscous timescale for the disc with ^ v = 0.1 at r p = 
1 AU: 




We consider the pressure effect of the gas disc and obtain 
the turnover time alo ng the horseshoe orbit in front of the 
planet at r v = 1 AU (|Baruteau fc Massed l2008h : 



Tturn — 




The viscous timescale is nearly equal to the turnover 
time in the disc with £ v = 0.1. The effect of the viscosity 
on the corotation torque then starts to be important. Note 
that the cooling timescale by the radiation is the same with 
the viscous timescale in the disc. 

The torque density near the orbit of the planet 
in the region 2 is small er than that in the disc B. 
iMorohoshi fc Tanakal ((2003) used a local shearing box cal- 
culation and studied the gravitational interactions between 



a planet and an optically thin disc, taking into account of 
energy dissipation by radiation. They found that the oval 
shape of the density contour profile near a planet is tilted 
with respect to the direction toward the central star. The 
asymmetry of the density structure increases the one-side 
torque. They showed that the density contour is less in- 
clined with the decreasing opacity. The adiabatic disc can 
be considered to have much larger opacity than the region 2, 
yiel ding the larger torque dens ity near the planet. 

IPaardekooper et al.l (|201 lh have derived the torque for- 
mula for the Lindblad torque and the corotation torque ex- 
erted on a planet by a disc with viscosity and thermal dif- 
fusion. We find that the corotation and Lindblad torques 
obtained in our simulation agree with those calculated from 
their torque formula within 30 %. The total torque is given 
by the sum of the Lindblad torque and the corotation torque. 
The total torque becomes Y = 2.1 in the region 1 and it is 
slightly smaller than that in the disc A, F = 3.0. The den- 
sity enhancement is reduced as the vortensity smears by the 
viscosity. However, the viscosity is too small to remove the 
corotatio n torque. Planets accumulate at the ice line as sug- 
gested bv lHasegawa fc Pudritj (|201l| y 

We further examine the effect of the viscosity of a disc 
on the torque density exerted on the planet, using the large 
value of the viscosity, £ v = 1-0. We utilize the same surface 
density distribution and the temperature distribution given 
by Fig.Q] Figure|H]shows the radial distribution of the torque 
density exerted on the planet in the region 1 at t = 25t p . The 
corotation torque is lost because the viscosity is large enough 
to smear the density enhancement. In this case, the viscous 
timescale becomes much shorter than the turnover time (see 
Eqs. (|30[) and l)3ip L The total torque becomes negative (r = 
—2.4), leading to the inward migration of the planet in the 
region 1. 

The density contours around the planet are shown to 
make clear the effect of the viscosity on the torque den- 
sity. Figure [9] shows the contours of the density fluctuations, 
E(i = 25tp)/Ei n i — 1.0, around the planet (a) in the disc A, 

(b) in the region 1 of the disc with (£ gr ,£v) = (10,0.1), and 

(c) in the region 1 of the disc with (£ gr ,£v) = (1.0,1.0). 
The enhancement of the density in the horseshoe orbit in 
the panel (b) is similar to that in the panel (a). On the 
other hand, the density in the horseshoe region is consid- 
erably reduced in the panel (c) because of the larger vis- 
cosity. The corotation torque decreases because the viscous 
timescale is much shorter than the turnover time in the case 
of £ v = 1.0. Moreover, the inner and outer spiral density 
waves are damped by the viscosity as well. This results in 
the small Lindblad torque. 

3.2.2 Dependence of the total torque on the opacity and 
the viscosity of a disc 

The corotation torque decreases in an adiabatic disc 
because the entropy of the gas in the horseshoe or- 
bit of the planet tends to b ecome uniform after the 
synodic period of the planet llBaruteau fc Massed 120081 ; 
IPaardekooper fc Papaloizoij|2008l . l2009af ). This is called the 
saturation of the corotation torque and happens in an adi- 
abatic disc. The corotation torque does not saturate in a 
viscous disc because the vortensity is transferred from the 
horseshoe region to the outer region of the disc by viscosity 
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IPaardekooper fc Papaloizoull2008l , l2009ar i. Figure [10] shows 
the time evolutions of the total torque acting on the planet 
in the adiabatic disc (the disc A) and in the optically thick 
accretion disc with (£ gr ,£v) = (10,0.1) (the region 1), which 
are represented by the solid and dotted curves, respectively. 
The total torque increases with time at the beginning of the 
simulation and reaches a steady state around t = 10i p . The 
total torque decreases in the adiabatic disc by the satura- 
tion, on the other hand, the total torque remains the same in 
the optically thick accretion disc. We do not have to worry 
about the saturation because we consider the planet in a 
disc with the viscosity. 

We examine the effect of the boundary of the regions 
of the disc on the total torque by changing the position 
of the planet. Figure [TTJ shows the total torques acting on 
the planets at the different positions in the optically thick 
accretion discs. The filled circles and triangles correspond to 
the total torques on the planet in the discs with (£ gr , £ v ) = 
(1.0, 1.0) and (10, 0.1), respectively. The total torque on 
the planet in the region 1 is nearly independent of f p . The 
effect of the boundary of the regions on the total torque 
is very small. The total torque in the region 2 decreases 
a little with f p . The analytical formula for the normalized 
corotation torque (|Paardekooper et al.ll201lh depends on the 
scale height and the viscous coefficient, and the corotation 
torque decreases with increasing f p in the region 2. This 
agrees with our numerical results. 

We perform a number of numerical simulations of grav- 
itational interactions between the planet and the optically 
thick accretion discs. Figure [12] shows the magnitude and 
sign of the total torques exerted on the planet in the re- 
gion 1 of the discs with various parameters of £ gr and 
£ v . The planet is located at 1 AU when fi2 = 2.4, 1.4, 
and 1.2, and at 0.7 AU when fi2 = 0.8, respectively. The 
open and filled marks denote the positive and negative to- 
tal torques acting on the planets, respectively. We obtain 
the largest magnitude of the total torques, |r| max , from the 
simulations with the same boundary position of fi2. The 
total torques are plotted as circles, squares, and triangles 
when |f| > 0.5|f | max , 0.5|f| max > |f| > 0.1|f| max , and 
0.1|F| max > |F|, respectively. The total torque increases with 
decreasing £ gr and £ v and its sign changes from negative to 
positive. 

We find the small decrease in the corotation torque 
when the viscous timescale is nearly equal to the turnover 
time as shown in Fig. 5. The corotation torque on the 
planet decreases with the viscosity of the disc and vanishes 
when the viscous timescale becomes much shorter than the 
turnover time as shown in Fig. 8. We find that the positive 
corotation torque cancels the negative Lindblad torque to 
have the zero total torque when r v i sc ~ 0.5rt ur n. We use the 
temperature distribution of the region 1 of the disc to find 
the viscous timescale: 




Using T v i sc = 0.5r tm . n , we obtain 



C gr = 3.8 x 10- 2 f^ ^ e v ~ 17/5 . (33) 

The dashed and dot curves in Fig. [12] are drawn using 
Eq. (|33[) . We consider the two positions of the planet, 



f p = 1.0 and 0.7. Substituting f p = 1.0 and 0.7 into Eq. ([33)) , 
we draw the dashed and dotted curves, respectively. The 
whole region is separated into the two regions by the curves. 
The positive total torque and the negative total torque can 
be found in one of the regions. We conclude that the coro- 
tation torque is moderately damped by the viscosity when 
T v i sc ~ 0.5T tU rn and cancels the negative Lindblad torque. 
The opacity of the disc increases with £ gr . Energy dissi- 
pates ineffectively in the disc with large £ gr , increasing the 
mid-plane temperature and the scale height. The viscous 
timescale of the disc with large £ gr decreases because the 
width of the horseshoe region decreases with the increasing 
scale height (see Eq. (|27|l 'l. In the disc with large £ gr , the vis- 
cosity quickly damps the corotation torque and the negative 
Lindblad torque becomes dominant. Moreover, it is valuable 
to mention the effect of the surface density at 1 AU, So, on 
Eq. (|33[) . In this study, we set Eo to be constant. Since the 
accretion rate and the optical depth are proportional to the 
surface density, r v i sc and r turn also increases with Eo. Hence, 
the coefficient of £ v on the right hand side of Eq. (|33[) changes 
by E . 

Additionally, one can find from Eq. (|33[) that the planet 
mass has a large effect on the threshold at which the total 
torque becomes zero. We examine the effect of the planet 
mass on the total torque. Figure [131 shows the total torque 
exerted on the planet with 3 Earth masses and that with 
7.5 Earth masses in the panels (a) and (b), respectively. The 
planet is located at f p = 1.0 in both cases. The half width 
of the horseshoe region increases with the planet mass. It 
takes longer time to damp the density enhancement in the 
horseshoe region of the planet with large mass. The larger 
viscosity is required to reduce the corotation torque acting 
on the massive planet. 



4 SUMMARY 

We have studied the Type I migration of a planet in an 
optically thick accretion disc. The gravitational interactions 
between a planet and disc gas excite spiral density waves 
inside and outside of the orbit of the planet in the disc. The 
waves attract the planet gravitationally and exert torques on 
the planet. The negative torque by the outer density wave is 
a little larger than the positive torque by the inner density 
wave because the outer wave is closer to the planet than the 
inner wave due to the negative pressure gradient. The sum 
of the torques by the density waves (the Lindblad torque) 
become s negative, leading to the inward migration of the 
planet (|Wardlll997n . 

The corotation torque is very important to planetary 
migration because it might be a ble to halt the inward mi- 
gration or reverse its direction (|Baruteau fc Massed 120081 ; 
iPaardekooper fc Papaloiz ou 2008). In a non-barotropic disc, 
it is found that the entropy related non-linear corota- 
tion torque plays an essential role. The non-linear coro- 
tation torque comes from the density enhancement due 
to the generated vort e nsity at the outgoing se paratrices 
(|Masset fc Casolill20ldl ; IPaardekooper et alJ[2oToh . The en- 
tropy related corotation torque is proportional to the radial 
gradient of the entropy of a disc. A disc with the steep en- 
tropy gradient yie l ds th e large positive corotation torque. 
lYamada fc Inaba I <|201lf > showed that the positive corota- 
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tion torque cancels the negative Lindblad torque when the 
power-law index of the entropy distribution of the disc be- 
comes A = —0.4, that is the critical power-law index of the 
entropy distribution. 

The radiation from the central star cannot reach the 
mid-plane of an optically thick disc. The temperature of the 
mid-plane of an optically thick accretion disc is determined 
by the energy balance between the viscous heating and the 
radiative cooling. The rate of the cooling by the radiation is 
sensitive to the opacity of the disc. Micron size dust particles 
provide the sources of the opacity. The opacity changes with 
the temperature through the sublimation of ice. We utilize 
the opacity law frequently used in other researches and find 
that the disc is divided into three regions: from the region 1 
with the highest temperature to the region 3 with the low- 
est temperature. The boundary positions of the regions are 
dependent on the viscosity and the opacity of the disc. The 
boundary between the regions 1 and 2 is located between 
0.5 AU and 2.5 AU. 

Each region has the different power-law indexes of the 
temperature and entropy distributions. The power-law in- 
dexes of the entropy distribution in the regions 1 and 3 are 
lower than the critical power-law index of —0.4, while the 
power-law index in the region 2 is higher. We found that 
the total torque exerted on the planet by the adiabatic disc 
becomes positive in the regions 1 and 3 and becomes nega- 
tive in the region 2, as expected. This means that the planet 
moves outward in the regions 1 and 3, while it moves inward 
toward the central star in the region 2. Planets might accu- 
mulate at the boundary between the regions 1 and 2 in the 
adiabatic disc. 

Dissipative processes chan ge the m ag nitude 



of the corotation torque (Kiev & Cridal 20081: 



Paardekooper fe Papaloizoull2008l . l2009al : lYamada fe Inaba I 
201ll ). In a viscous disc, the vortensity is transferred 



radially as shown in Fig. [S] The viscosity is responsible 
for the decrease in the corotation torque. The radiation 
has an impact on the evolution of the entropy as well and 
keeps the entropy gradient within the horseshoe region. The 
corotation torque depends on these diss ipative processes 
due to the viscosity and the radiation (|Masset fe Casohl 
l20ld ). The positive total torque exerted on the planet 
by an adiabatic disc might become negative once some 
dissipative processes are included. We include the viscosity 
and radiation into the adiabatic disc and calculate the total 
torque exerted on the planet in an optically thick accretion 
disc. We focus on the planet in the regions 1 and 2. The 
magnitudes of the opacity and the viscosity are expressed 
with the two parameters, £ gr and £ v , respectively. The 
opacity and the viscosity increase with £ gr and £ v - The large 
mass accretion rate of the disc is found in the disc with large 
£ v . The temperature of the disc increases when the energy 
of the disc dissipates ineffectively (large £ gr ), resulting in 
the narrow horseshoe region. The viscosity decreases the 
enhancement of the density in the horseshoe region. The 
corotation torque gets smaller when the viscosity is larger. 
The total torque exerted on the planet by the optically 
thick accretion disc depends on the two parameters. 

We have made a number of numerical simulations of 
gravitational interactions between the planet and the op- 
tically thick accretion disc with various parameter sets 
of £ gr and £ v - The total torque always becomes negative 



in the region 2, leading to the inward migration of the 
planet. The total torque becomes positive in the region 1, 
if the effect of the dissipation is small. The dissipative 
processes work effectively and the total torque becomes 
zero when the timescale for the viscosity is half of the 
turnover time of the planet in the horseshoe orbit, r v i sc = 
0.5rt U rn- This equation is written with the parameters 
as 5v = 0.4^ 1/17 (Af p /5Af E ) 18 / 17 (Eo/100g/cm 2 ) 7 / 17 e gr 5/17 , 
considering Eo as a parameter. In the optically thick accre- 
tion disc with f; gr ~ 1.0 and the surface density of 100 g/cm 2 
at 1 AU, the accretion rate of the disc is required to be 
smaller than 2.1 x lO _8 M0/yr for the planet to move out- 
ward. Our study suggests that planets with 5 Earth masses 
might accumulate and drive further growth of the planets 
at the boundary between the regions 1 and 2 in the opti- 
cally thick accretion disc. The small accretion rate of gas is 
required for small planets to move outward. 
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The analytical formula suggests that the Lindblad torque 
increases with a decrease in a and agrees with numerical 
results. The corotation torque becomes a little larger in the 
disc with a = 0.5, compared with that in the disc with a = 
1.0. This increase in the corotation torque is also suggested 
by the analytical formula for the corotation torque formula. 
We can derive the relation between £ gr and £ v when the total 
torque becomes zero, using the same procedure explained in 
the section 3.2. Using TVisc = 0.5T turn , we have 

fc = «« io-,-"" (&) "" (jj^jj) "" «.-"".(A3) 



APPENDIX A: EFFECT OF a ON THE TOTAL 
TORQUE ON A PLANET 

We examine the dependence of the total torque on the 
power-law index, a, of the surface density distribution of 
the disc. Figure \A1\ shows the time evolutions of the to- 
tal torque exerted on the planet by the disc with a — 1.0 
and the disc with a = 0.5. Both discs have parameters of 
(CgnCv) = (10,0.1). The solid and dashed curves represent 
the total torques acting on the planet by the discs with 
a = 1.0 and a — 0.5, respectively. The total torques in- 
crease at the beginning of the simulations and eventually 
reach steady states in both discs. The total torque by the 
disc with a = 0.5 becomes approximately 1.5 times as large 
as that by the disc with a = 1.0. 

The torque densities by the disc with a = 1.0 and by 
the disc with a = 0.5 are plotted as the solid and dashed 
curves, respectively, in Fig. \A2\ The inner and outer Lind- 
blad torques increase and decrease with a decrease in a, 
respectively, leading to the smaller Lindblad torque. The 
analytical formula for the Lindblad torque, Pl, is given by 
|Paardekooper et al]|201ll ) 

7 f L = -2.5- 1.7ft + 0.1a. (Al) 

The temperature distribution of the region 1 of the optically 
thick accretion disc has the power-law index, j3i = (3 + a)/3. 
Hence, we have for the Lindblad torque in the region 1: 



7 r L = -4.2 - 0.47a. 



(A2) 
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Figure 1. (a) the surface density distribution and (b) the temperature distribution of the disc with a = 1.0. The temperature distribution 
of the disc is determined by the balance between the viscous heating and the radiative cooling. The temperature distribution is described 
by the power-law distributions with three different gradients. The power-law index of the temperature distribution strongly depends on 
the opacity of the disc. The dependence of the opacity on the temperature changes due to the evaporation of ice. The power-law indexes 
of the temperature distribution change at 1.4 AU and 2.5 AU. 
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Figure 2. The boundary position between the regions 1 and 2, fi2, as a function of £ gr and £ v . The solid, dashed, dot-dashed, and 
dotted lines correspond to fia = 2.4, 1.4, 1.2, and 0.8, respectively. We use disc models with 20 parameter sets (§gr,£v) marked by the 
filled circles. 
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Figure 3. The radial distributions of the torque density acting on the planet with 5 Earth masses by the adiabatic disc at f = 25t p . 
The power-law indexes of the temperature distributions of the disc are /3 = 4/3 and /3 = A/11 in the panels (a) and (b), respectively. 
The planets are located at 1 AU and 2 AU in the panels (a) and (b), respectively. The horizontal axis is normalized by r p . The strong 
corotation torque is found in the panel (a), leading to the outward migration of the planet. 
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Figure 4. The half width of the horseshoe region as a function of the planet mass. The horizontal axis is normalized by the Earth mass. 
The optically thick disc with the accretion rate of £ v = 0.1 is used. The numerical half width of the horseshoe region (filled circles) is 
evaluated at t = 25i p . T he planet is located at f D = 1.0. The curve is drawn using the analytical expression for the half width of the 
horseshoe region given by [Paardckoop er et alj J2010T) . 
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Figure 5. Same as Fig. 3, but the optically thick accretion disc with (§ gr ,Cv) = (10, 0.1). The normalized total torque acting on the 
planet is 2.1 in the panel (a), while it becomes —1.5 in the panel (b). The planet migrates outward in the region 1 and inward in the 
region 2. Planets might accumulate at the boundary of the regions. 
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Figure 6. Contour of the vortensity in (a) the inviscid adiabatic disc (disc A) and (b) the accretion disc with the non- vanishing viscosity 
((£v, £gr)=(0.1,10)) at t = 10tp. The initial distribution of vortensity is subtracted from the vortensity for better contrast. The horizontal 
and vertical axes correspond to r/r p and 9/n, respectively. The planet with 5 Earth masses is located at 1 AU. 
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Figure 7. The corotation torques acting on the planet with 5 Earth masses given by the analytical expression of lMasset &: Casol i (2010) 
(solid line) and the numerical results (filled circles). The horizontal and vertical lines represent the viscous strength factor, £ v , and the 
normalized corotation torque, r c , respectively. In deriving the numerical corotation torque, we used Eq. H27t for the half width of the 
horseshoe region. The planet is located at 1 AU. 




Figure 8. The radial distribution of the torque density acting on the planet with 5 Earth masses by the region 1 of the disc with 
(£gr,£v) = (1.0, 1.0). The planet is located at 1 AU. The total torque acting on the planet becomes negative, V = —2.4. Planets move 
inward toward the central star even in the region 1. 
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Figure 9. The contour plots of the surface density enhancement T,(t = 25t p )/Si n i — 1.0 around the planet with 5 Earth masses (a) in an 
adiabatic disc (disc A), (b) in the region 1 of the disc with (§gr,£v) = (10, 0.1), and (c) in the region 1 of the disc (£ gr ,£v) = (1.0, 1.0), 
respectively. The horizontal axis is normalized by r p and the vertical axis is divided by 7r. The planet is located at (r, 8) = (r p , 7r). The 
density contour in the panel (b) is similar to that in the panel (a) because the viscosity is too small to decrease the density enhancement. 
On the other hand, the density enhancement is greatly reduced by the viscosity in the panel (c). 
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Figure 10. The time evolution of the total torque acting on the planet with 5 Earth masses in the adiabatic disc (solid curve) and in 
the region 1 of the disc with (£gr,£v) = (10, 0.1) (dashed curve). In both cases, the planet is located at 1 AU. The total torque decreases 
around t = 25t p in the adiabatic disc due to the saturation of the corotation torque. On the other hand, the positive total torque is 
maintained in the disc with (£ g r,£v) = (10, 0.1). 
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Figure 11. The total torques exerted on the planet by the disc with (£ g r,§v) = (1.0, 1.0) (filled circle) and (10, 0.1) (filled triangle) as 
a function of f p . The boundary position between the regions 1 and 2 is shown by the dotted line. In the region 1, the total torque seems 
to be nearly independent of f p . In the region 2, the total torque decreases a little with increasing f p . 
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Figure 12. The diagram of the total torques exerted on the planet with 5 Earth masses by the disc with various parameter sets of 
(£gr,£v). The open and filled marks correspond to the positive and negative total torques acting on the planets, respectively. We find 
the largest magnitude of the total torques, |r| max , from the simulations with the same boundary position of fia- The total torques are 
plotted as circles, squares, and triangles when T 0.5|r| max , 0.5|r| max > |r| > o.i|r| max; and 0.1|r| m ax > |T|, respectively. The dotted 
and dashed curves are drawn using T v j sc — O.STturn when 7'p — 0.7 and 1.0, respectively. It is noted that the curves are truncated because 
f\2 is smaller than f p . The sign of the total torque changes when r v i sc ~ 0.5rt ur n- The planet with 5 Earth masses moves outward in 
the region 1 of the optically thick accretion disc with the surface density of 100 g/cm 2 at 1 AU when the accretion rate is smaller than 
2.1 X 10 _8 MQ/yr. The planets are then expected to accumulate at the boundary of the regions 1 and 2 in the optically thick disc. 
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Figure 13. Same as Fig, 1121 but the planet with (a) 3 Earth masses and (b) 7.5 Earth masses. The open and filled marks correspond 
to the positive and negative total torques, respectively. The dotted curve is drawn using r v ; sc = 0.5Tt urn when f p = 1.0. The planets 
with 3 and 7.5 Earth masses moves outward in the region 1 of the optically thick accretion disc with the surface density of 100 g/cm 2 
at 1 AU when the accretion rate becomes smaller than 1.2 X 10~ 8 and 3.2 X 10~ 8 Mg/yr, respectively. 
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Figure Al. The time evolution of the total torque exerted on the planet with 5 Earth masses by the disc with (£ gr , §v) = (10, 0.1). The 
solid and dashed curves correspond to the total torques acting on the planet by the discs with a = 1.0 and a = 0.5, respectively. The 
planet is located at 1 AU. The total torque by the disc with a = 0.5 is larger than that by the disc with a = 1.0. 
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Figure A2. The torque densities exerted on the planet by the discs with a = 1.0 (the solid curve) and a = 0.5 (the dashed curve). 
Both discs have the parameters of (£ gr ,£v) = (10, 0.1). The planet is located at 1 AU. 



